1 Les techniques bootstrap
1.1 Introduction
Les techniques bootstrap reposent sur une une idée similaire. Ici, nous allons échantillonner, avec remise, de manière répétée parmi les individus de notre échantillon et nous allons approximer la distribution d’échantillonnage d’une statistique donnée. L’idée clé ici est que, même si la solution du modèle de régression est déterministe, les données elles-mêmes sont supposées être échantillonnées de manière aléatoire, et donc toutes les estimations que nous faisons dans un modèle de régression sont aléatoires. Si les données changent, alors les estimations changent aussi. Le bootstrap nous donne une appréciation non-paramétrique de la distribution de ces estimations. Une fois encore, l’avantage de cette méthode est que nous pouvons construire des intervalles de confiance pour, par exemple, la pente de la droite de régression, sans avoir à supposer que les résidus sont normalement distribués.
## Loading required package: Stat2DataDans cet exemple, nous examinons les prix des Porsche. Nous voulons estimer le prix (Price) en fonction du kilométrage (Mileage). Construisons une représentation graphique de notre échantillon de 30 voitures :
Pour créer un échantillon bootstrap, nous sélectionnons des lignes de la base de données uniformément au hasard, mais avec remplacement.
xyplot(Price ~ Mileage, data=mosaic::resample(PorschePrice), pch=19, cex=2, alpha=0.5, type=c("p","r"), lwd=5)
- Notez les différences entre ce graphique et le précédent. Que remarquez-vous ?
- Exécutez le code précédent plusieurs fois. Qu’est-ce qui change ? Qu’est-ce qui reste inchangé ?
1.2 Interactions et réflexions sur le bootstrap
Le package manipulate nous permet de créer des
échantillons bootstrap à la volée.
bslope <- function (formula, data, n) {
# Data d'origine
# Mise en oeuvre du bootstrap
bootstrap <- mosaic::do(n) * coef(lm(formula, data=mosaic::resample(data)))
xyplot(formula, data=data,
panel = function (x, y, ...) {
panel.xyplot(x,y, pch=19, cex=2, alpha=0.5)
fm <- lm(formula, data=data)
panel.abline(fm, col="red", lwd=5)
# Ajout des pentes bootstrap
for (i in 1:n) {
panel.abline(t(bootstrap[i,]), -0.5, col="lightgray", lwd=0.3)
}
panel.text(75, 80, paste("mean intercept\n", round(mean(~Intercept, data=bootstrap), 6)), cex=0.75)
panel.text(75, 75, paste("sd intercept\n", round(sd(~Intercept, data=bootstrap), 6)), cex=0.75)
panel.text(75, 70, paste("mean slope\n", round(mean(~Mileage, data=bootstrap), 6)), cex=0.75)
panel.text(75, 65, paste("sd slope\n", round(sd(~Mileage, data=bootstrap), 6)), cex=0.75)
}
)
}
# bslope(Price ~ Mileage, data=PorschePrice, 100)Expérimentez avec l’exemple suivant jusqu’à ce que vous ayez une idée de ce qui se passe.
if(!require(manipulate)){install.packages("manipulate");library(manipulate)}
manipulate(
bslope(Price ~ Mileage, data=PorschePrice, n),
n = slider(2, 500, initial=100)
)
1.3 Distributions bootstrap et intervalles de confiance
L’un des avantages du bootstrap est qu’il nous permet de construire une distribution d’échantillonnage pour le coefficient de pente qui ne dépend pas du respect des conditions d’utilisation du modèle de régression linéaire.
Les intervalles de confiance originaux de notre modèle RLS (régression linéaire simple) dépendent de la véracité de ces conditions.
fm <- lm(Price ~ Mileage, data=PorschePrice)
confint(fm) # Astuce, si vous voulez changer le niveau de confiance, consultez l'aide## 2.5 % 97.5 %
## (Intercept) 66.2360186 75.9448869
## Mileage -0.7054401 -0.4733618Maintenant, créons une distribution bootstrap pour les coefficients de la régression.
Nous utilions seulement 100 échantillons mais il faudrait en utiliser plus !
bootstrap <- mosaic::do(100) * coef(lm(Price ~ Mileage, data=mosaic::resample(PorschePrice)))
p2 <- densityplot(~Mileage, data=bootstrap)
ladd(panel.abline(v=coef(fm)["Mileage"], col="red", lwd=3), plot=p2)
La distribution bootstrap sera toujours centrée sur la valeur de nos données réelles, mais elle nous indique d’autres valeurs probables pour le coefficient (essentiellement, l’erreur d’échantillonnage). Une façon de quantifier cette variabilité est de créer un intervalle de confiance.
Dans la suite, nous présentons trois méthodes pour construire des intervalles de confiance à partir du bootstrap. La méthode la plus robuste est la méthode n°2.
1.3.1 Méthode n°1
La première méthode suppose que la distribution bootstrap est normale. (C’est une hypothèse étrange. Si elle était normale, nous utiliserions simplement une distribution normale pour approximer la distribution d’échantillonnage). Dans ce cas, nous pourrions utiliser l’écart-type de la distribution bootstrap pour construire un intervalle de confiance pour le coefficient de pente.
## [1] -0.6840458 -0.49475611.3.2 Méthode n°2
La seconde méthode n’exige pas que la distribution bootstrap soit normale, mais fonctionne mieux lorsqu’elle est grossièrement symétrique. Dans ce cas, nous utilisons simplement les percentiles de la distribution bootstrap pour construire les intervalles de confiance. Cette méthode est la plus logique dans la plupart des cas.
## 2.5% 97.5%
## -0.6758628 -0.49053011.3.3 Méthode n°3
Si la distribution bootstrap est biaisée, nous devons adapter la deuxième méthode pour qu’elle fonctionne encore. En effet, notre estimation peut déjà différer du véritable paramètre de population.
Dans le cas particulier de notre exemple, la distribution bootstrap pour le coefficient de pente n’est pas biaisée.
qs <- qdata(~Mileage, p = c(0.025, 0.975), data=bootstrap)
coef(fm)["Mileage"] - (qs - coef(fm)["Mileage"])## 2.5% 97.5%
## -0.5029390 -0.6882718- Créez un échantillon bootstrap d’au moins 2000 et construisez les trois intervalles de confiance mentionnés ci-dessus. En quoi diffèrent-ils de l’intervalle de confiance typique ?
- Construisez des intervalles de confiance bootstrap pour les
coefficients des modèles
GPA∼SATVetGPA∼SATM.
2 Tests de randomisation (Permutation)
2.1 Introduction
Rappelons que l’inférence pour la régression linéaire repose sur l’hypothèse que les trois conditions suivantes sont remplies :
- Linéarité
- Variance constante
- Normalité des résidus
Nous savons comment évaluer si ces conditions sont remplies, et vous avez vraisemblablement appris quelques techniques pour les corriger lorsqu’elles ne le sont pas (par exemple, les transformations). Aujourd’hui, nous apprendrons quelques techniques basées sur la randomisation pour faire des inférences non paramétriques. Ces inférences ne reposent pas sur des hypothèses strictes concernant la distribution des résidus.
Dans cet exemple, nous essayons d’expliquer la moyenne générale d’un étudiant de l’université en fonction de son score à l’épreuve orale SAT.
Commençons par explorer le jeu de données.
library(DT)
datatable(FirstYearGPA, rownames = FALSE, filter="top", options = list(pageLength = 5, scrollX=T) )Examinons d’abord la relation entre ces deux variables à l’aide d’un graphique.
if(!require("lattice")){install.packages("lattice");library(lattice)}
lattice::xyplot(GPA ~ SATV, data=FirstYearGPA, pch=19, cex=2, alpha=0.5, type=c("p","r"), lwd=5)
Il est aussi possible de créer un graphique interactif.
library(ggplot2)
#library(plotly)
p <- FirstYearGPA %>%
ggplot( aes(SATV, GPA, size = HSGPA, color=Male)) +
geom_point() +
stat_smooth(aes(x = SATV, y = GPA), method = "lm", formula = y ~ x, se = F, color="gold") +
stat_smooth(method = "loess", formula = y ~ x) +
scale_x_log10() +
theme_bw()
plotly::ggplotly(p)## Warning: The following aesthetics were dropped during statistical transformation: size.
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
## the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
## variable into a factor?## Warning: The following aesthetics were dropped during statistical transformation: size
## and colour.
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
## the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
## variable into a factor?Il semble y avoir une certaine association linéaire entre ces variables, mais elle n’est pas particulièrement forte. Nous pouvons quantifier cette relation en utilisant le coefficient de corrélation.
## [1] 0.3043114Dans ce cas, la valeur du coefficient de corrélation est d’environ 0,3, ce qui n’est pas grand, mais semble être significativement différent de 0. Rappelez-vous le \(t\)-test de corrélation que vous avez appris dans le passé :
\[ t=r\sqrt{\frac{n-2}{1−r^2}}, \quad n−2,\ \textrm{ddl}. \]
##
## Pearson's product-moment correlation
##
## data: x and y
## t = 4.706, df = 217, p-value = 4.499e-06
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.1789575 0.4199430
## sample estimates:
## cor
## 0.3043114Ce test confirme que la corrélation entre la moyenne générale à l’universté et la SATV n’est probablement pas nulle, mais la validité de ce test exige que certaines hypothèses soient respectées (plus précisément, la normalité bivariée du couple dont la corrélation est testée). Supposons que dans ce cas, les hypothèses ne soient pas satisfaites.
Pouvons-nous encore être sûrs que la corrélation est non nulle ?
2.2 Permuter les données
Si la moyenne générale et la SATV étaient réellement corrélées, alors il existe une relation réelle liant la ième valeur de la moyenne générale à la ième valeur de la SATV Dans ce cas, il ne serait pas logique de lier la ième valeur de la moyenne générale à une autre valeur de la SATV. Mais si la corrélation entre ces deux variables était en fait nulle, alors la façon dont nous faisons correspondre les entrées dans les variables n’aurait aucune importance !
L’idée de base du test de permutation est de mélanger plusieurs fois la correspondance entre les deux variables (c’est-à-dire l’échantillon sans remplacement), et d’examiner la distribution du coefficient de corrélation qui en résulte. Si la valeur réelle du coefficient de corrélation est une réalisation rare de cette distribution, alors nous aurons la preuve que l’hypothèse nulle d’absence de corrélation (corrélation nulle) pourrait ne pas être vraie.
- Exécutez le code suivant plusieurs fois. Vous pouvez ainsi vous faire une idée de la façon dont la ligne de régression peut changer, en fonction de la permutation des données.
- Croyez-vous toujours que la corrélation de la pente est non nulle ?
La procédure pour le test de randomisation est simple. Il suffit de mélanger la variable de réponse et de calculer le coefficient de corrélation résultant avec la variable explicative. Mais nous le faisons plusieurs fois, jusqu’à ce que nous ayons une idée de la distribution de ce coefficient de corrélation. Nous examinons alors où se situe la corrélation observée dans cette distribution.
# Répétons 5000 fois l'opération
rtest <- mosaic::do(5000) * cor(GPA ~ shuffle(SATV), data=FirstYearGPA)
p1 <- densityplot(~cor, data=rtest, xlim=c(-0.4,0.4), xlab="Correlation Coefficient")
ladd(panel.abline(v=cor.actual, col="red", lwd=3), plot=p1)
Bien sûr, nous pouvons aussi trouver explicitement où se situe dans la
distribution la corrélation observée.
## [1] 0Le paramètre lower.tail indique à R
s’il doit calculer la probabilité cumulée à droite ou à gauche de la
réalsation de la statistique du test. Dans ce cas, nous voulons calculer
la probabilité cumulée à droite de la ligne, donc nous fixons
lower.tail=FALSE.
Enfin, nous pouvons trouver un intervalle de confiance non paramétrique de 95% pour le coefficient de corrélation. L’interprétation qui peut en être faite ici est que si notre valeur réelle de corrélation tombait dans cet intervalle, nous ne la considérerions pas comme statistiquement significative au seuil de 5%.
## 2.5% 97.5%
## -0.1330114 0.1324548Comparez cet intervalle de confiance à celui renvoyé par le test de corrélation ci-dessus. Pourquoi sont-ils différents ?
Effectuez la procédure ci-dessus pour tester la corrélation entre la GPA et la SATM.
Effectuez la procédure ci-dessus pour tester la corrélation entre la GPA et la HSGPA.
2.3 Pour aller plus loin
- Article de Tim Hesterberg publié en 2015 dans The American Statistician et intitulé “What Teachers Should Know About the Bootstrap: Resampling in the Undergraduate Statistics Curriculum”.
- Article fondateur de Brad Efron